    function [tabMC,tabchi,tabtchi,theored]=conflevel(Mspecred,specred,nsim,nw,Ax,npr,dt,rho,fr)
i=1;
% estimate of the MC (Markov chains) confidence levels.
n85=round(nsim*85/100);
n90=round(nsim*90/100);
n95=round(nsim*95/100);
n99=round(nsim*99/100);

tabvecsort=[];
% sort the power values of specred for each frequency
for i=1:npr
    vec=specred(i,:);
    tabvecsort(i,:)=sort(vec);
end  

tabMC=[];
for i=1:npr
    tabMC(i,1)=tabvecsort(i,n85);
    tabMC(i,2)=tabvecsort(i,n90);
    tabMC(i,3)=tabvecsort(i,n95);
    tabMC(i,4)=tabvecsort(i,n99);
end 

% Estimate of the Chi-square confidence levels for the mean spectrum of the 
% nsim red noise signals.

tabchi=[];

% The degrees of freedom for a mtm analysis is equal to 2*number of tapers
% (Tompson, 1982).
nw2=2*(2*nw-1);
facchi85=chi2invPMTK(0.85,nw2)/nw2;
facchi90=chi2invPMTK(0.90,nw2)/nw2;
facchi95=chi2invPMTK(0.95,nw2)/nw2;
facchi99=chi2invPMTK(0.99,nw2)/nw2;

tabchi(:,1)=Mspecred*facchi85;
tabchi(:,2)=Mspecred*facchi90;
tabchi(:,3)=Mspecred*facchi95;
tabchi(:,4)=Mspecred*facchi99;

% estimate of the theoretical red noise spectrum
fnyq=1/(2*dt);
theored=[];
i=1;
for i=1:npr
   theored(i)=(1-rho^2)/(1-(2*rho*cos(pi*fr(i)/fnyq))+rho^2);
end;
% normalisation of the spectrum
theoredun=theored(1);
theored(1)=0;
Art=sum(theored)/npr;
theored(1)=theoredun;
theored=theored*(Ax/Art);

% Estimate of the Chi-square confidence levels theoretical red noise
% spectrum.

tabtchi=[];

tabtchi(:,1)=theored*facchi85;
tabtchi(:,2)=theored*facchi90;
tabtchi(:,3)=theored*facchi95;
tabtchi(:,4)=theored*facchi99;
end